library(memisc)                         #spss.portable.file

setwd("code/")
source("insurance_source_score_1_prep.R")
setwd("code/")
source("insurance_source_score_3_functions.R")


hni130 <- as.data.set(spss.portable.file(paste0("data/hni130.por")))
hni130$vulnerable <- factor(
    hni130$Q20,
    levels=c("Feel vulnerable to high medical bills","Feel well-protected by your health insurance plan")
)
## hni130$PSRAID <- hni130$psraid
hni130$NUMBER <- 130

hni138 <- as.data.set(spss.portable.file(paste0("data/p2015hni138.por")))
hni138$vulnerable <- factor(
    hni138$Q41,
    levels=c("Feel vulnerable to high medical bills","Feel well-protected by your health insurance plan")
    )
## hni138$PSRAID <- hni138$psraid
hni138$NUMBER <- 138

hni149 <- as.data.set(spss.portable.file(paste0("data/hni149.por")))
hni149$vulnerable <- factor(
    hni149$Q18,
    levels=c("Feel vulnerable to high medical bills","Feel well-protected by your health insurance plan")
    )
## hni149$PSRAID <- hni149$psraid
hni149$NUMBER <- 149

sdat.under65.2_add <- merge(
    sdat.under65.2,
    rbind(
        data.frame(hni130[,c("PSRAID","NUMBER","vulnerable")]),
        data.frame(hni138[,c("PSRAID","NUMBER","vulnerable")]),
        data.frame(hni149[,c("PSRAID","NUMBER","vulnerable")])
        ),
    by=c("PSRAID","NUMBER"),
    all.x=T
)


getVar <- function(file_number, var, thelevels=c("Your health insurance premiums","The deductible you pay before insurance kicks in","Co-pays for doctor visits and prescription drugs","Some other health care cost","Or is paying for health care and health insurance not a financial burden for you?", "All equally (VOL.)"),varname="burden") {
    hni <- as.data.set(spss.portable.file(paste0("data/",file_number,".por")))
    hni[,varname] <- factor(data.frame(hni[,var])[,var])
    ## hni$burden <- factor(data.frame(hni[,var])[,var], levels=thelevels)
    ## hni$PSRAID <- hni$psraid
    hni$NUMBER <- gsub("hni", "", file_number)
    return(hni)
}
hni100 <- getVar("hni100", "Q23")
hni110 <- getVar("hni110", "Q15")
hni142 <- getVar("hni142", "Q27")
hni151 <- getVar("hni151", "Q7")
hni154 <- getVar("hni154", "Q4")

sdat.under65.2_add <- merge(
    sdat.under65.2_add,
    rbind(
        data.frame(hni100[,c("PSRAID","NUMBER","burden")]),
        data.frame(hni110[,c("PSRAID","NUMBER","burden")]),
        data.frame(hni142[,c("PSRAID","NUMBER","burden")]),
        data.frame(hni151[,c("PSRAID","NUMBER","burden")]),
        data.frame(hni154[,c("PSRAID","NUMBER","burden")])
        ),
    by=c("PSRAID","NUMBER"),
    all.x=T
)

lms_uninsured_acafavor <- plot_slopes_return_lmer_list(
    outcome="FAVOR>=3",
    propensity="COVERED == 0",
  the_subset=subset(
    sdat.under65.2_add,
    INCOME > 40
  ),
  plot_title="Uninsured Score\nrespondents under 65, income over 40K",
  pdf_name="scratch.pdf",
    color=gray(0.2),
  no_pdf = TRUE
)

df <- lms_uninsured_acafavor[["out_data"]]



df$burden[df$burden %in% c("(DO NOT READ) Don't know","(DO NOT READ) Refused","(DO NOT READ) Don't know/Refused")] <- NA

## table(subset(df, !is.na(burden))$DATE)
## table(df$burden)

mod_copays <- (
    lm(
        as.integer(burden%in%c(
            "Your doctor visits",
            "Your prescription drugs",
            "Co-pays for doctor visits and prescription drugs"
        )) ~
            PSCORE3,
        data = subset(df, !is.na(burden)),
        weights = WEIGHT
    )
)

mod_premiums <- (
    lm(
        as.integer(burden%in%c(
            "Your health insurance premiums"
        )) ~
            PSCORE3,
        data = subset(df, !is.na(burden)),
        weights = WEIGHT
    )
)

mod_deductibles <- (
    lm(
        as.integer(burden%in%c(
            "The deductible you pay before insurance kicks in"
        )) ~
            PSCORE3,
        data = subset(df, !is.na(burden)),
        weights = WEIGHT
    )
)

mod_no_burden <- (
    lm(
        as.integer((burden%in%c(
            "Or is paying for health care and health insurance not a financial burden for you?"
        ))) ~
            PSCORE3,
        data = subset(df, !is.na(burden)),
        weights = WEIGHT
    )
)

stargazer(
    mod_premiums,
    mod_copays,
    mod_deductibles,
    mod_no_burden,
    type="latex",
    digits=2,
    column.labels=c("Premiums","Co-pays","Deductibles","No burden"),
    star.cutoffs = c(0.05, 0.01, 0.001),
    out=paste0(
        the_prefix, "/figs/",
        "tableA23_insurance_score_burden_models.tex"
    )
)


mod_vulnerable <- (
    lm(
        I(vulnerable == "Feel vulnerable to high medical bills") ~
            PSCORE3,
        data = subset(df, !is.na(vulnerable)),
        weights = WEIGHT
    )
)

stargazer(
    mod_vulnerable,
    type="latex",
    digits=2,
    star.cutoffs = c(0.05, 0.01, 0.001),
    out=paste0(
        the_prefix, "/figs/",
        "tableA24_insurance_score_vulnerable_model.tex"
    )
)
